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ABSTRACT: Numerical studies of dislocation pair correlations have played a cen- 
tral role in deriving a continuum theory from the equations of motion of 2D dislocation 
systems in a mathematically rigorous way. As part of an effort to extend this theory 
into the full 3D dislocation problem, 3D dislocation pair correlations were studied 
with discrete dislocation dynamics simulation. As a first approximation, dislocations 
were modeled as uncharged curves in space (their Burgers vectors were disregarded). 
An inverse square decay with distance was found to describe the numerically obtained 
pair correlations of the studied curve system. 
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1 INTRODUCTION 

Statistical mechanics studies of discrete 2D systems of straight, parallel edge dislo- 
cations have led Groma and co-workers to a rigorously derived prototype continuum 
model for bridging the micro- and mesoscales of crystal plasticity. Numerical calcu- 
lation of dislocation pair correlation functions played a crucial part in the construction 
of the theory 0101: i) The observed short range nature of the correlations (which is 
nontrivial due to long range elastic dislocation interactions) is a prerequisite for a local 
theory to exist, ii) Pair correlation functions of homogeneous systems are explicitly 
present in the theory in the flow stress term as well as in a non-local diffusion-like 
term which enables the theory to describe size effects. 

In current efforts to generalize this theory to the full 3D plasticity problem (31 0], 
dislocation pair correlations are expected to play a similar role as in the prototype 
2D theory [5]. As there exist no theoretical predictions on the functional form of 3D 
dislocation correlations, numerical simulation is the obvious means to provide the nec- 
essary input for the continuum theory. A pioneering 3D discrete dislocation dynamics 
(DDD) simulation study [6] considered the deformation of a symmetrically oriented 
large bcc crystal with periodic boundary conditions and analyzed the simulated dislo- 
cation network as a system of randomly distributed curves (in mathematical terms, a 
stochastic fiber process), disregarding the Burgers vector information. Unfortunately, 
the single shot approach of this study did not provide enough statistics to get smooth 
graphs as results and the authors did not connect their numerical results to the crucial 
theoretical questions mentioned above. 

In the present paper, the results of a large number of statistically equivalent DDD sim- 
ulations are presented. As in 0, the obtained dislocation configurations are analyzed 
as a system of randomly distributed uncharged curves in space and the radial decay of 
the pair correlation function of this curve system is studied. 



2 SIMULATIONS 



The 3D dislocation configurations analyzed below were generated with the DDD code 
described in 0, [SI]. This code is optimized for the small scale plasticity of finite, 
cuboid-shape fee single crystals and does not handle periodic boundary conditions. 
In the model, dislocations are discretized into connected straight segments which are 
allowed to glide. Junction formation and dislocation annihilation upon contact are 
also included. Cross slip was turned off in the present simulations to mimic low 
temperature conditions. 

In the following we present the simulated uniaxial tensile deformation of (0.8 um) 3 Al 
cubes up to an applied strain e a ppiied = 0.67% with an applied strain rate 5000 s -1 . The 
strain rate was chosen just below the quasistatic limit, i.e. above which the measured 
stress-strain curves started to show enhanced hardening due to inertia effects. To 
improve angular statistics, the crystals were oriented for [010] symmetric multiple 
slip. The side faces of the specimens were traction free. At the top and bottom faces 
displacement boundary conditions were prescribed in the tensile direction whereas 
the in plane components were traction free, too. The simulations were started with 16 
randomly positioned and oriented Frank-Read sources of length 0.22 um in each slip 
system, resulting in an initial dislocation density pinitiai = 8.3 x 10 13 mT 2 . 
As mentioned in the Introduction, for the present analysis the simulated dislocation 
configurations were simplified by considering all dislocations and junctions as sim- 
ple curves in space (their Burgers vectors were disregarded). We then calculated the 
radial pair correlation function g(r) of these curve systems (see eq. fl} for its defini- 
tion) which was complicated by three practical problems, i) To get smooth curves as 
results, we needed to improve the statistics by averaging over an ensemble of identi- 
cal simulations, different only in their random initial configurations. In the following, 
an ensemble of 55 simulations is analyzed which proved large enough to smooth out 
both the average hardening behavior (see fig. HQl and the averaged pair correlation 
functions (see below), ii) As can be seen in fig. [2] the obtained dislocation configura- 
tions were inhomogeneous with dislocation depleted zones near the sample surfaces. 
The thickness of these zones, around 0.2 um, was found to be roughly equal to the 
mean dislocation-dislocation distance (see also fig.[T|for volume averaged val- 

ues of p). Despite the different boundary conditions at the top and side cube faces, 
no visible differences were detected between these two types of dislocation depleted 
zones. Therefore, in the following we restrict our analysis to the central (0.4 um) 3 
sub-volume which proved homogeneous to a good approximation (see fig. O. iii) 
The linear size of the studied sub-volume was only four times larger than the mean 
dislocation-dislocation distance 1 / yfp. This distorted each obtained pair correlation 
function g(r) already at distances r where it was still far from its asymptotic value. We 
corrected for this effect by dividing the measured pair correlation functions with the 
pair correlation function of a random point distribution in the (0.4 um) 3 sub-volume. 



3 DISCUSSION 

The radial pair correlation function of a random system of curves can be defined as 



||9t] where L(r) is the mean total line length in a sphere of radius r centered at a "typ- 
ical" curve point and p means the line density (total length in a unit volume). In this 
context, "typical" means that a random choice is made in a way that every point on 
the curves has the same chance to be chosen. The numerical method to compute g(r) 
directly followed its definition: i) The entire dislocation network (both mobile disloca- 
tions and junctions) was divided into line segments of length 0.5 nm. ii) For every line 
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Figure 1: Solid line: typical simulated stress-strain curve. Dashed line: ensemble 
averaged stress-strain curve. Dotted line: evolution of the ensemble and volume av- 
eraged dislocation density. 



segment the total curve length was calculated in a sphere of radius r. iii) The curve 
length values were averaged over all line segments and then substituted into eq. (Q]) as 
the numerical estimate for L(r). The results of such a calculation at e a ppiied = 0.67% 
are displayed in fig. [3] (solid line). (It was found that g(r) depends only slightly on 
eappiied and it seems to saturate with increasing e a ppiied- Further simulations to confirm 
this up to eappiied = 1.3% are in progress). As can be seen in fig. [3] g(r) diverges 
as r~ 2 as r — ► and converges to 1 as r — ► oo. The main goal of this study is to 
establish the rate of this convergence. Before this can be done, however, g(r) needs 
further analysis for two reasons, i) The contiguity of dislocation lines causes strong 
correlations at small r values. These need to be separated from g(r) as only the cor- 
relations between different dislocation lines are interesting for theory 0]. ii) In the 
simulations, dislocation lines are discretized into linear segments of length 13-66nm 
which might also distort the numerically obtained g(r) data. In the following, we 
analyze g(r) along these problems. 

The connected nature of dislocation lines causes two types of strong correlations at 
small r: i) self correlations of smooth dislocation curves g se \f(r) and ii) correlations 
between different dislocations near their intersections (e.g. arms going out of junc- 
tions) <7ints(?")- Fi rst we analytically calculate the r — ► asymptotes of these two 
contributions utilizing eq. £[), then separate them from the numerical g(r) data at 
small r values to estimate their impact. 

As r — > 0, <7 se ]f (r) is equivalent to the pair correlation function of randomly distributed 
straight lines with a linear density p. We only need that part of the correlation function 
where a line correlates with itself. From eq. (Q~|) it readily follows that 



(the same result was also derived in [9]). 

The r — > asymptote of gi nts (r) can be modeled with "corners" homogeneously dis- 
tributed in space. In a "corner" two arms of length I meet at an angle (p and the arms 
have a linear density p. For the r — > asymptote of <?i nts , only that part of the cor- 
relation function is interesting where an arm correlates with the arm attached to it. 
Straightforward calculation from eq. (fl} yields 
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Figure 2: Evolution of the ensemble averaged density profile, averaged over the tensile 
and one of the perpendicular directions. 

Note that in the case of simulated dislocation configurations, averaged I and <p values 
appear in eq. (f3]) but this does not affect the r _1 character of gi ats (r) at small r. 
The next step is assessing the relative contributions of g se \t(r) and <?mts(?~) to the nu- 
merically calculated pair correlation function g(r). To this end, we separated g(r) 
into three terms utilizing the fact that the simulation code stores dislocation segments 
organized into loops. The first term, g se \f(r), was computed as the correlation function 
of segment pairs residing on the same loop. The second one, flints Mi was calculated 
as the pair correlation of segment pairs on different loops which touch each other. The 
last term, g pa ir(r), was defined simply as fl(r) — fl se if(r)— flints (r). Fig.|3]displays these 
three contributions and the total pair correlation function g(r) at e app ij e d = 0.67%. Nu- 
merical fits of the r~ 2 asymptote of g se if(r) and the r" 1 asymptote of flints(f) at r — > 
are also displayed. Both fit well for r < 10 nm, consistently with the minimum dis- 
cretization segment length 13nm, and diverge from the numerical curves at larger 
distances. Note that fl P air(?") shows no singularity at r — > 0, i.e. all connectivity related 
pair correlation terms were successfully separated. Note also that fl pa ir(r) < 1 (is an- 
ticorrelated) and that, to a smaller extent, even g(r) < 1 for r > 0.1 um (not visible 
on the graph). The latter is clearly a finite size effect caused by the independence of 
Jo°° g( r ) r2 dr from fl(r), that g(r) = 1 is a valid pair correlation function, and that 
g{r) ^> 1 as r — > 0. The authors took the liberty to manually compensate for this 
effect by adding a small constant c -C 1 to g(r) when analyzing its decay to 1. 
We finally study the decay of g(r) towards the uncorrected value 1 as r — ► oo. Fig. [4] 
depicts g(r) — 1 + c for e app iied = and e ap piied = 0.67%. With c = (not shown in 
the figure), a g(r) — 1 oc r~ 2 decay can be seen up to r « 20 nm for e app ii e d = and 
r as 30 nm for e app ii e d = 0.67% before g(r) submerges 1. However, the authors think 
such an anticorrelated g(r) unlikely. Instead, by carefully adjusting c, an amazingly 
good fit to r~ 2 can be achieved: at e app iied = with c = 0.25 up to r = 50 nm and at 
eappiied = 0.67% with c = 0.08 up to r = 100 nm (1.5 times the discretization segment 
length 66 nm; see fig.|U). It is clear from fig.[3]that at e app iied = 0.67% this r~ 2 decay 
is not plainly a result of the same-loop correlation <? se if(r) (cf. eq. (f2|) as <? S eif(?~) is 1.5 
orders of magnitude smaller than g(r) at r = 100 nm. It has also nothing to do with 
individual straight segments as they are shorter than 66 nm. Therefore, the observed 




Figure 3: Decomposition of the pair correlation function g(r) at 0.67% applied strain 
(see text for details). The r — ► asymptotes of <?seif('") oc r~ 2 and gi nt s(r) oc r" 1 , and 
the r — > oo asymptotes of g(r) — * 1 and g p - d i r (r) — > 1 are also indicated. 



r 2 decay of g(r) is a collective effect at larger distances. This is notable from the 
point of view that the upper limit of this decay in our simulations, r = 100 nm, is still 
in the range of the mean dislocation-dislocation distance 1 / y/p (see fig. Q}. 




Figure 4: Evolution of the r — > oo decay of the pair correlation function g(r). 



In summary the results indicate a gradual extension of an r 2 correlated zone from 
small r values towards larger ones as the deformation proceeds. A new simulation 
effort is ongoing to reach e a ppiied = 1-3% for further confirmation. 
In a final note, a pair correlation function oc r~ 2 should be proportional to p^ 1 for 
dimensional reasons. This is consistent with the numerical results in fig. [4] (see also 
fig-fflfor P values). 



4 CONCLUSIONS 



The correlation properties of 3D many-dislocation systems were studied with discrete 
dislocation dynamics simulation in symmetrically oriented fee crystals deformed in 
uniaxial tension. As a first approximation, dislocation configurations were analyzed 
as uncharged curve systems in space and the evolution of the corresponding radial 
pair correlation function was studied. The results indicate the gradual appearance of a 
pair correlation function ( o _1 r -2 from small r towards larger ones as the deformation 
proceeds. Interestingly, this decay is identical to that of 2D systems of infinite, parallel 
edge dislocations relaxed in multiple slip [10]. Moreover, substituting the numerical 
result /9 _1 /|r — r'| 2 into eqs. (10-11) in J3|] as the double angular integral of d pml ~ 
suggests that the dimensionless constants denoted by a and D(9) in [5] have at most 
a logarithmic divergence with system size. Whether this divergence is cancelled by 
the combined angular dependence of g and the dislocation pair interaction r, enabling 
finite a and D(9) values and thus a local 3D dislocation continuum theory, is the 
objective of subsequent numerical work. 
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